capture log close
clear all
set maxvar 10000
set more off
pause on
# delimit;

local crime "";
local output "";


/* 
RD estimates for "exactage" dataset
*/

local exp "prworaexact";

local fs_vars  "age18_yes n_age18unfav n_onssi1922 n_onssi1922_exact n_onssi20 " ;
local hh_vars "n_earnhh9702 n_earnhh9712 n_emphhpost9702 n_emphhpost9712 n_inchh9702 n_inchh9712" ;
local covariates "n_male n_firstage n_fam_singmom n_fam_nopar n_dibmdr_mie n_dibmdr_mine n_dibmdr_miss n_diag1_mental n_diag1_nerv n_diag1_sense n_diag1_infec n_diag1_endo n_diag1_none n_diag1_cong n_diag1_musc n_diag1_resp n_diag1_blood n_diag1_neo n_earnhhpre9095" ;
local new_vars "ct_*" ;
local robust_vars "n_crime2char_18_38 n_crime2chart_18_38 n_crime2genr_18_38 n_crime2genrt_18_38  n_crime2nogn_18_38 n_crime2nognt_18_38 n_crime2ahinc_18_38" ;

local controls_all "n_male i.n_firstage i.n_diag1cat n_dibmdr_mie n_dibmdr_mine n_dibmdr_miss n_fam_singmom n_fam_nopar i.n_earnhhbins i.state" ;
local controls_no  "n_dibmdr_mie n_dibmdr_mine n_dibmdr_miss" ; 
local controls_qje1 "" ; 
local controls_qje2 "n_male i.n_firstage i.n_diag1cat n_fam_singmom n_fam_nopar i.n_earnhhbins i.state" ;

local bwidth_list "264" ;
local cov_list   "all no" ; 
local demo_main "all male female" ;
local demo_oth "firstlt6 first7to13 firstgt13 singno twopar parhigh parlow othmen nonment intel phys nonphys white black nonw 
				ue97high ue97low uegrhigh uegrlow ueallhigh uealllow" ;
local demo_hat "q25 q2575 q75" ;

foreach demo in `demo_main' { ;
	local vars_`demo' "n_crime* c_crime* n_1crime* n_2crime* n_3crime* n_4crime* n_yr* avg_* ct_*" ; 
	local othvars_`demo' "`covariates' `fs_vars' `hh_vars'" ;
} ;
foreach demo in `demo_oth' { ;
	local vars_`demo' "n_crime*" ;
	local othvars_`demo' "`fs_vars'" ;
} ;
foreach demo in `demo_hat' { ;
	local vars_`demo' "n_crime2char_18_38 n_crime2chart_18_38 n_crime2genr_18_38 n_crime2genrt_18_38  n_crime2nogn_18_38 n_crime2nognt_18_38 n_crime2ahinc_18_38" ;
	local othvars_`demo' "" ;
} ;

local poly_list "1(1)1" ;

use "`crime'\prep2\crime_prwora_prepforreg_exactage.dta", clear ;

/* Create predicted removal probability */
xi: reg n_age18unfav 
	`controls_all' i.dibdig
	if n_post == 1 ;
predict hat_removal ;

keep if abs(dobdd_run) <= `bwidth_list' ;

/* Subsamples */
gen demo_all = 1 ;
gen byte demo_othment	= (n_diag1_mental == 1 & (substr(dibdig,1,3) != "317" & substr(dibdig,1,3) != "318" & substr(dibdig,1,3) != "319" )) ;
gen byte demo_nonment 	= (demo_othment == 0) ;
gen byte demo_intel 	= (n_diag1_mental == 1 & (substr(dibdig,1,3) == "317" | substr(dibdig,1,3) == "318" | substr(dibdig,1,3) == "319" )) ;
gen byte demo_phys 		= (n_diag1_mental != 1) ;
gen byte demo_nonphys 	= (n_diag1_mental == 1) ;
gen byte demo_twopar 	= (n_fam_singmom != 1 & n_fam_nopar != 1) ;
gen byte demo_singno 	= (n_fam_singmom == 1) | (n_fam_nopar ==1) ;
gen byte demo_nopar  	= (n_fam_nopar == 1) ;
gen byte demo_male 	 	= (n_male ==1) ;
gen byte demo_female 	= (n_male !=1) ;
gen byte demo_white		= (race == "W") ;
gen byte demo_black 	= (race == "B") ;
gen byte demo_nonw		= (race != "W") ;

summ hat_removal, det ;
gen byte demo_q25 	= (hat_removal <  r(p25)) ;
gen byte demo_q2575 = (hat_removal >= r(p25) & hat_removal < r(p75)) ;
gen byte demo_q75 	= (hat_removal >= r(p75) & hat_removal != .) ;

egen cnty_ue97 = rowmean(cnty_urate_1997-cnty_urate_1998) ;
egen cnty_uegr = rowmean(cnty_urate_2007-cnty_urate_2012) ;
egen cnty_ueall = rowmean(cnty_urate_1997-cnty_urate_2017) ;
rename cnty_crimerate cnty_crime ;
rename cnty_share_blue cnty_dem ;
rename cnty_pop2000 cnty_pop ;

foreach cnty in crime dem pop ue97 uegr ueall { ;
	summ cnty_`cnty', det ;
	gen byte demo_`cnty'high = (cnty_`cnty' > r(p50) & cnty_`cnty' != .) ;
	gen byte demo_`cnty'low = (cnty_`cnty' <= r(p50) & cnty_`cnty' != .) ;
	replace demo_`cnty'high = . if cnty_`cnty' == . ;
	replace demo_`cnty'low = . if cnty_`cnty' == . ;
	tab demo_`cnty'high, m ;
	tab demo_`cnty'low, m ;
} ;

summ n_earnhhpre9095, det ;
replace n_earnhhpre9095 = 0 if n_earnhhpre9095 == . ;
gen byte demo_parhigh = (n_earnhhpre9095 > r(p50) & n_earnhhpre9095 != .) ;
gen byte demo_parlow = (n_earnhhpre9095 <= r(p50) & n_earnhhpre9095 != .) ;
replace demo_parhigh = . if n_earnhhpre9095 == . ;
replace demo_parlow = . if n_earnhhpre9095 == . ;
tab demo_parhigh, m ;
tab demo_parlow, m ;

gen byte demo_record = (n_crime2any_pre18 == 1) ;

gen byte demo_firstlt6 		= (n_firstage <=6) ;
gen byte demo_first7to13 	= (n_firstage > 6 & n_firstage <=13) ;
gen byte demo_firstgt13 	= (n_firstage > 13 & n_firstage < .) ;

/* RD variables */
gen dobdd_run2=dobdd_run*dobdd_run;
gen dobdd_run3=dobdd_run2*dobdd_run;
gen dobdd_run4=dobdd_run3*dobdd_run;

gen postXdobdd=n_post*dobdd_run;
gen postXdobdd2=n_post*dobdd_run2;
gen postXdobdd3=n_post*dobdd_run3;
gen postXdobdd4=n_post*dobdd_run4;

local poly_1 "dobdd_run postXdobdd";
local poly_2 "dobdd_run dobdd_run2 postXdobdd postXdobdd2";
local poly_3 "dobdd_run dobdd_run2 dobdd_run3 postXdobdd postXdobdd2 postXdobdd3";
local poly_4 "dobdd_run dobdd_run2 dobdd_run3 dobdd_run4 postXdobdd postXdobdd2 postXdobdd3 postXdobdd4";

foreach demo in all { ;

foreach bwidth in `bwidth_list' {;

	forval poly = `poly_list' {;
		
		foreach cov in `cov_list' { ;
		
				
				/* RD regressions */
				
				reg n_post n_post `poly_`poly'' 
					if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1 
					, robust;
					outreg2 using "`output'\\`exp'_bw`bwidth'_crime_poly`poly'`cov'cov_rd_`demo'.xls", e(N r2 rmse) replace;
				
				foreach var of varlist 
					`vars_`demo''
					{;
					xi: reg `var' n_post `poly_`poly''
						`controls_`cov''
						if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1
						, robust; 
					outreg2 using "`output'\\`exp'_bw`bwidth'_crime_poly`poly'`cov'cov_rd_`demo'.xls", e(N r2 rmse) append;
					} ;	
				
				foreach var of varlist 
					`othvars_`demo''
					{ ;
			
					xi: reg `var' n_post `poly_`poly''
						`controls_`cov''
						if abs(dobdd_run)<=`bwidth' & keep_1997_2017_2c == 1 & demo_`demo' == 1
						, robust; 
					outreg2 using "`output'\\`exp'_bw`bwidth'_crime_poly`poly'`cov'cov_rd_`demo'.xls", e(N r2 rmse) append;
				} ; /* end var */		
			

		
		/* IV estimates */
				
				ivreg2 n_post (n_age18unfav=n_post) `poly_`poly''
					if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1
					, robust;
					outreg2 using  "`output'\\`exp'_bw`bwidth'_crime_poly`poly'`cov'cov_iv_`demo'.xls", e(N r2 rmse) replace;
				
				/* IV using unfavorable age 18 */
				foreach var of varlist 
					`vars_`demo'' 
					{;
					
					qui summ `var' 
						if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1 ;
					if r(sd) != 0 { ;
						xi: ivreg2 `var' (n_age18unfav=n_post) `poly_`poly''
							`controls_`cov''
							if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1
							, robust;
						outreg2 using  "`output'\\`exp'_bw`bwidth'_crime_poly`poly'`cov'cov_iv_`demo'.xls", e(N r2 rmse) append;
						} ;
					};
				
				if "`demo'" == "all" | "`demo'" == "male" | "`demo'" == "female" { ;
					foreach var of varlist 
						`hh_vars' 
						{ ;
					
						qui summ `var' 
							if abs(dobdd_run)<=`bwidth' & keep_1997_2017_2c == 1  & demo_`demo' == 1 ;
						if r(sd) != 0 { ;
							xi: ivreg2 `var' (n_age18unfav=n_post) `poly_`poly''
								`controls_`cov''
								if abs(dobdd_run)<=`bwidth' & keep_1997_2017_2c == 1  & demo_`demo' == 1
								, robust;
							outreg2 using  "`output'\\`exp'_bw`bwidth'_crime_poly`poly'`cov'cov_iv_`demo'.xls", e(N r2 rmse) append;
						} ;
					} ; /* end var */
				} ; /* end if */
				
			} ; /* end cov */
				
		} ; /* end poly */
	} ; /* end bwidth*/
} ;	/* end demo */


/* Look at thresholds */
local thresh_list "char genr" ;
local thresh_max = 20 ;

local ct_list "" ;
foreach bwidth in `bwidth_list' {;				
foreach demo in all { ;
forval poly = `poly_list' {;
foreach cov in `cov_list' { ;
		
	foreach stub in `thresh_list' {;
		
		foreach thresh in exact gt { ;

			forval i = 1(1)`thresh_max' { ;

			qui summ ct_`thresh'_crime2`stub'`i'
				if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1 ;
			if r(sd) != 0 { ;
				local ct_list "`ct_list' `thresh'`stub'`i'" ;
				xi: statsby _b _se, saving("`output'\statsby\\`exp'_ct`thresh'`stub'`i'_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta", replace):
					 ivreg2 ct_`thresh'_crime2`stub'`i' (n_age18unfav=n_post) `poly_`poly''
					`controls_`cov''
					if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1
					, robust;
				} ;	/* end if */
			
			} ; /* end i */
		} ; /* end thresh */
	} ; /* end stub */
	
} ; /* end cov */
} ; /* end poly */
} ; /* end demo */
} ; /* end bwidth */

foreach bwidth in `bwidth_list' {;				
foreach demo in all { ;
forval poly = `poly_list' {;
foreach cov in `cov_list' { ;
				
			use 			"`output'\statsby\\`exp'_ctexactchar1_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta", clear ;
			gen var = "exactchar1" ;
		
			foreach var in `ct_list' { ;
				
				append using "`output'\statsby\\`exp'_ct`var'_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta" ;
				replace var = "`var'" if var == "" ;
					
				} ;
			} ;
			drop if _n == 1 ;
			
			save 			"`output'\\`exp'_ctcompiled_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta", replace ;
		
} ;
} ;
} ;



/* Any/number of charges by certain date */
local cumul_varlist " c_*char* c_*genr* c_combpay*" ; /*  */
local cumul_varlist_full "" ;

foreach bwidth in `bwidth_list' {;				
foreach demo in all { ;
forval poly = `poly_list' {;
foreach cov in `cov_list' { ;

	matrix means_`bwidth'_`demo'_`poly'_`cov'=J(200,1,.) ;
	local count = 0 ;
		
	foreach var of varlist `cumul_varlist' {;
	
		summ `var' if n_post == 0 & dobmo_run >= -6 & dobmo_run < 0 & demo_`demo'== 1 ;
		
		if r(mean) > 0 { ;
		
			local count = `count' + 1 ;
			local cumul_varlist_full "`cumul_varlist_full' `var'" ;
			
			/* Save control mean */
			matrix means_`bwidth'_`demo'_`poly'_`cov'[`count',1]=r(mean) ;
			
			/* Save regression estimate */
			xi: statsby _b _se, saving("`output'\statsby\\`exp'_cumul_`var'_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta", replace):
				 ivreg2 `var' (n_age18unfav=n_post) `poly_`poly''
				`controls_`cov''
				if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1
				, robust;
		} ;

	} ; /* end var */
	
} ; /* end cov */
} ; /* end poly */
} ; /* end demo */
} ; /* end bwidth */

foreach bwidth in `bwidth_list' {;				
foreach demo in all { ;
forval poly = `poly_list' {;
foreach cov in `cov_list' { ;
clear ;

	svmat means_`bwidth'_`demo'_`poly'_`cov' ; 
	rename means_`bwidth'_`demo'_`poly'_`cov'1 mean ;
	gen var = "" ;
	local new_count = 0 ;
	foreach var in `cumul_varlist_full' { ;
		local new_count = `new_count' + 1 ;
		replace var = "`var'" if _n == `new_count' ;
	} ;
	drop if mean == . ;
	save "`output'\statsby\\`exp'_cumul_means.dta", replace ;

} ;
} ;
} ;
} ;

foreach bwidth in `bwidth_list' {;				
foreach demo in all { ;
forval poly = `poly_list' {;
foreach cov in `cov_list' { ;
				
			use 			"`output'\statsby\\`exp'_cumul_c_crime2char18_18_38_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta", clear ;
			gen var = "c_crime2char18_18_38" ;
		
			foreach var in `cumul_varlist_full' { ;
				
				append using "`output'\statsby\\`exp'_cumul_`var'_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta" ;
				replace var = "`var'" if var == "" ;	
			} ;
			drop if _n == 1 ;
			keep var _b_n_age18unfav  _b_dobdd_run  _b_postXdobdd  _b_cons
					 _se_n_age18unfav _se_dobdd_run _se_postXdobdd _se_cons ;
			
			merge 1:1 var using "`output'\statsby\\`exp'_cumul_means.dta" ; 
			save 			"`output'\\`exp'_cumul_compiled_`bwidth'_poly`poly'`cov'cov_iv_`demo'.dta", replace ;
		
} ;
} ;
} ;
} ;




/* Robustness to different bandwidths */
foreach var of varlist `robust_vars' {;
				
	foreach demo in all { ;

		forval poly = `poly_list' {;
			
			foreach cov in `cov_list' { ;
	
				forval bwidth = 25(25)400 { ;
				
					xi: statsby _b _se, saving("`output'\statsby\\`exp'_robust`bwidth'_`var'_poly`poly'`cov'cov_rd_`demo'.dta", replace):
						reg `var' n_post `poly_`poly''
							`controls_`cov''
							if abs(dobdd_run)<=`bwidth' & demo_`demo' == 1
							, robust; 
					} ;	
			
			} ;
		} ;
	} ;
} ;

foreach var of varlist `robust_vars' {;
				
	foreach demo in all { ;

		forval poly = `poly_list' {;
			
			foreach cov in `cov_list' { ;
			
				use 			"`output'\statsby\\`exp'_robust25_`var'_poly`poly'`cov'cov_rd_`demo'.dta", clear ;
				gen bwidth = 25 ;
				
				forval bwidth = 25(25)400 { ;
					append using "`output'\statsby\\`exp'_robust`bwidth'_`var'_poly`poly'`cov'cov_rd_`demo'.dta" ;
					replace bwidth = `bwidth' if bwidth == . ;
				} ;
				drop if _n == 1 ;
				
				order bwidth _b_n_post _se_n_post ;
				save 			"`output'\\`exp'_robust_`var'_poly`poly'`cov'cov_rd_`demo'.dta", replace ;
				
			} ;
		} ;
	} ;
} ;


capture log close;
